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ABSTRACT 

We discuss the nature of the various modes of pulsation of superfluid neutron 
stars using comparatively simple Newtonian models and the Cowling approx- 
imation. The matter in these stars is described in terms of a two-fluid model, 
where one fluid is the neutron superfluid, which is believed to exist in the 
core and inner crust of mature neutron stars, and the other fluid represents a 
conglomerate of all other constituents (crust nuclei, protons, electrons, etc.). 
In our model, we incorporate the non-dissipative interaction known as the 
entrainment effect, whereby the momentum of one constituent (e.g. the neu- 
trons) carries along part of the mass of the other constituent. We show that 
there is no independent set of pulsating g-modes in a non-rotating superfluid 
neutron star core, even though the linearized superfluid equations contain a 
well-defined (and real- valued) analog to the so-called Brunt- Vaisala frequency. 
Instead, what we find are two sets of spheroidal perturbations whose nature 
is predominately acoustic. In addition, an analysis of the zero-frequency sub- 
space (i.e. the space of time-independent perturbations) reveals two sets of 
degenerate spherioidal perturbations, which we interpret to be the missing 
g-modes, and two sets of toroidal perturbations. We anticipate that the de- 
generacy of all these zero-frequency modes will be broken by the Coriolis force 
in the case of rotating stars. To illustrate this we consider the toroidal pul- 
sation modes of a slowly rotating superfluid star. This analysis shows that 
the superfluid equations support a new class of r-modes in addition to those 
familiar from, for example, geophysical fluid dynamics. Finally, the role of 
the entrainment effect on the superfluid mode frequencies is shown explicitly 
via solutions to dispersion relations that follow from a "local" analysis of the 
linearized superfluid equations. 
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1 INTRODUCTION 

In the forty or so years since Migdal's ( |1959| ) initial work, three important areas of research 
on superfluidity in neutron stars have emerged: (i) Nuclear physics studies of the pairing gap 
energies as functions of mass density (for recent reviews, see (|Lombardo 1999[ [Lombardo & 
Schulze 2000|) ); (ii) observations and subsequent modeling of glitches (see, for instance, ( [Sauls 
1989| ; |Lyne 19931) and references therein) and neutron star cooling flTsuruta 1994j ; [Heiselberg 



8z Hjorth- Jensen 2000|) ; and (iii) studies of the dynamics of superfluid neutron star cores in 
the Newtonian QEpstein 198l| |Mendell fc Lindblom 19911 ; |Mendell 19911 ; [Lindblom fc Mendell 
19941 ; [Lindblom fc Mendell 19951 ; [Lee 1995| ; [Prix 1999| ; |Lindblom fc Mendell 200"0|) and general 



relativistic regimes ( Carter 1989) ; |Comer fc Langlois 1994| ; Carter fc Langlois 1995a| ; Carter 
Langlois 1995b] ; Carter fc Langlois 1998| ; |Langlois et al 199S| ; Comer et al 1999] [Andersson 
fc Comer 200 1| ). It is this third area that is the focus of the present discussion. We want to 
understand how the pulsation properties of a neutron star changes when it cools below the 
temperature (a few times 10 9 K) at which the bulk of the interior is expected to become 
superfluid. One of the underlying motivations for this work is the possibility that future 
gravitational-wave detectors may be sensitive enough to observe pulsating neutron stars, 
following (say) a glitch or a starquake ([Andersson fc Kokkotas 1996| ; |Andersson fc Kokkotas 
1998 ; Kokkotas et al 2001|) . An analysis of such observed data will provide a probe of the 



neutron star interior and should, at least in principle, allow us to infer the details of the 
supranuclear equation of state. In particular, one would hope to be able to establish beyond 
any doubt that the core of a neutron star contains a superfluid and perhaps constrain some 
of the relevant parameters for (say) entrainment. Obviously we can only hope to achieve 
this goal if the pulsation properties of a superfluid star differ significantly from those of an 
"ordinary fluid" neutron star. 

The equations that describe ordinary (cold, inviscid etcetera) fluid neutron stars consist 
of a single mass continuity equation and a single Euler equation that determines the fluid 
velocity. But even this "simple" system yields an amazing diversity of modes of pulsation 
(see McDermott et al ( |1988| ) for a fine discussion of many classes of modes), which includes 
the spheroidal f-, p-, and g-modes and the toroidal r-modes. In reality cold neutron stars are 
significantly more complicated than what is implied by the ordinary fluid model. The current 
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On the dynamics of superfluid neutron star cores 3 

understanding is that the outer regions of a mature neutron star consist of crust nuclei 
and an electron gas, which (beyond neutron drip) are everywhere permeated by superfluid 
neutrons. In the core, nuclei have dissolved leaving behind neutrons in a superfluid state, 
superconducting protons, and an ultra-relativistic gas of degenerate electrons. In addition, 
more exotic particles (e.g. hyperons) may be present. Deep in the core various condensates 
(of kaons, pions etcetera) may form, and deconfined quarks may also play a crucial role. A 
priori, one might expect the pulsation spectrum of such a "real" neutron star to be very 
complex, consisting of the various modes that exist for the ordinary fluid case, plus additional 
modes that arise because of new fluid degrees of freedom due to different particle species 
moving (more or less) independently of each other. 

The most striking ways that (pure) superfluids differ from ordinary fluids is that su- 
perfluids have zero viscosity and are locally irrotational ( [Tilley fc Tilley 1986| ). The latter 
property, however, is compensated by the superfluid being threaded by an array of quantized 
vortices so that it can on macroscopic scales mimic the rotational behaviour of an ordinary 
fluid. The neutron superfluid in a rotating neutron star is believed to contain such threading. 
The vortices also lead to an effective viscosity through the so-called entrainment effect and a 
consequence of it known as mutual friction. In a mixture of the two superfluids Helium three 
and Helium four, it is known that a momentum induced in one of the constituents will cause 
some of the mass of the other to be carried along, or entrained flAndreev fc Bashkin 1975 ; 
Tilley fc Tilley 19861 ). The analog in neutron stars is an entrainment of some of the protons, 
say, by the neutrons ([Alpar et al 1984 ; Borumand et al 1996 ). Because of entrainment, the 
flow of neutrons around the neutron fluid vortices will also induce a flow in a fraction of the 
protons, leading to magnetic fields being formed around the vortices. But since the electrons 
are coupled to the protons on very short timescales ([Alpar et al 1984| ), some will track closely 
the entrained protons. Mutual friction is the dissipative scattering of these electrons off of 
the magnetic fields associated with the vortices. 

The comparatively simple model of neutron star superfluidity that will be used in this 
paper considers just two fluids, and the entrainment effect that acts between them. One 
fluid consists of the superfluid neutrons that exist in the inner crust and core, whereas the 
other is a conglomerate of all the charged particles (i.e. the nuclei and electrons in the crust, 
and the superconducting protons and electrons in the core) that will be loosely referred to 
as "protons." In comparison to the ordinary fluid case, the superfluid system of equations 
will consist of two mass continuity equations, and two Euler equations that determine the 
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neutron and proton fluid velocities. This model has been used in previous numerical studies 
of the linearized pulsations of superfluid neutron stars ( |Lindblom fc Mendell 1994] ; [Lee 19951 ; 
Corner et al 1999| ). However, it is important to point out that Mendell ( |1991| ) has discussed 



the more general case where electromagnetic effects are explicitly included, and where the 
electrons (and even muons) are free to move independently of the protons, a net result of 
which is to delimit the dynamical timescales for which the simplified two-fluid model is 
appropriate. 

It is also relevant to mention that even before the equations were analyzed in detail, 
Epstein ( |1988| ) (and then later Mendell (|1991|) ) argued, using a simple counting of the fluid 
degrees of freedom in a superfluid neutron star, that there should be a new class of modes 
that were later dubbed superfluid modes. They arise primarily because the neutrons and 
protons in the core are superfluid and no longer locked together in nuclei, thus leading to 
an increase in the number of fluid degrees of freedom. It is also impressive that — using as 
an analogy a system of coupled pendulums — Mendell ( |1991| ) argued that the distinguishing 
characteristic of the new modes should be the counter-motion of the neutrons with respect to 
the protons. That is, in terms of a projection along the radial direction, the superfluid modes 
should have the protons moving oppositely to the neutrons, unlike an ordinary fluid mode 
that would have the neutrons and protons moving more or less in "lock-step." This picture 
has been confirmed by numerical (and simplified analytical) studies in both the Newtonian 
( [Lindblom fc Mendell 19941 ; [Lee 1995| ) and general relativistic regimes flComer et al 1999| ) . 
Another feature, shown by Lindblom and Mendell ( 1994|) , is that the superfluid modes are 
driven by deviations from chemical equilibrium between the neutrons and protons. 

At first sight, the existence of the superfluid modes appears to confirm one's intuition 
that the new fluid degrees of freedom would double the number of pulsation mode-families. 
However, there is an open question in the literature regarding the g-modes in the superfluid 



case. In particular, Lee's ( |1995|) numerical results did not reveal any independent set of pul- 
sating g-modes, much less a doubling. This is perhaps strange, and certainly not consistent 
with the simple counting argument given above, since Reisenegger and Goldreich ( |1992| ) 
have shown convincingly that a composition gradient (for instance, a stable stratification in 
terms of the proton fraction) leads to the existence of propagating g-modes in non-rotating 
cold ordinary fluid neutron stars. Even more puzzling is the fact that Lee showed that the 
ordinary fluid definition of the so-called Brunt- Vaisala frequency leads to real values, which 
is the usual indication that pulsating g-modes should exist. His final conclusion was that 
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the superfluid modes are to be found among the ordinary fluid f- and p-modes, and that 
there are no propagating g-modes in the cores of superfluid neutron stars, even if there ex- 
ists a composition gradient. Here, we will use a "local" analysis of the linearized superfluid 
equations to show that superfluid modes are predominately acoustic and to confirm Lee's 
numerical results that superfluidity prevents g-modes from pulsating. 

Of course, if there are no distinct g-modes then one is still left with the question of 
where they have gone. To find the answer, we will analyze the zero-frequency subspace of 
solutions to the superfluid equations of motion, following recent studies by Lockitch and 
Friedman ( |1999|) and Lockitch et al ([2001 ) in the ordinary fluid case. This subspace is 
composed of the time-independent solutions to the linearized equations. We find that the 
zero-frequency subspace appropriate to the superfluid equations is spanned by solutions that 
separate into two distinct classes: (i) Non-trivial local scalar matter (e.g. number densities, 
pressure, etc) perturbations accompanied by zero velocity perturbations, and (ii) degenerate 
velocity perturbations with vanishing local matter perturbations. The first class of solutions 
simply takes static and spherically symmetric configurations and deforms them into other, 
nearby static and spherically symmetric configurations. The second class separates into 
independent sets of spheroidal and toroidal velocity perturbations. This is exacly analogous 
to the results for an ordinary fluid, but there is one important difference. We find two sets 
of decoupled perturbations for each of the spheroidal/toroidal velocity perturbations. 

We interpret these results in the following way: Because the second class of solutions 
have vanishing local matter variations, it is natural to take the spheroidal solutions to be 
the missing g-modes. The existence of two independent sets of such solutions would then 
be in agreement with the intuitive notion that the fluid degrees of freedom are doubled 
in the superfluid case. The presence of spheroidal modes in the zero-frequency subspace 
is analogous to the ordinary fluid case when there is no stratification in the background 
star ( |Lockitch fc Friedman 1999| ; [Lockitch et al 200 1|) . The doubling of modes in the class 
of toroidal solutions is also unique to the superfluid. Given these results, an analysis of 
rotating configurations becomes of prime importance. A comparison to the ordinary fluid 
case ( [Lockitch fc Friedman 19U9"| ; |Lockitch et al 20TJT )) suggests that rotation will break the 
degeneracy of the various zero-frequency and lead to inertial modes with true dynamics. Our 
results for the zero-frequency subspace suggests that, just like in a barotropic ordinary fluid 
star, the typical inertial mode will be a hybrid mixture of spheroidal and toroidal velocity 
components in the non-rotating limit. 
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6 N. Andersson and G. L. Comer 

We take the first steps towards the study of inertial modes in a superfTuid star by consider- 
ing the so-called r-modes (which correspond to purely toroidal velocities in the non-rotating 
limit). These modes are of particular relevance since it was recently discovered ([Andersson 
|1998| ; [Friedman fc Morsink 1998[ ) that they are generically unstable to gravitational radiation 
via the so-called CFS (Chandrasekhar-Friedman-Schutz) mechanism ( |(Jhandrasekhar 1970| ; 
Friedman fc Schutz 1978j ; [Friedman 1978| ). In principle, the CFS mechanism can apply to any 
mode in a rotating neutron star, but viscous damping of various kinds tends to kill the insta- 
bility before substantial amounts of gravitational radiation are emitted, except, apparently, 
for the r-modes ( [Lindblom et al 1998| ; |Andersson et al 1999| ; |Owen et al 1998] ; [Andersson] 
|fc Kokkotas 200l[) . Remarkably, even though mutual friction in superfluid neutron stars is 
extremely effective at killing the CFS instability for f- and p-modes, Lindblom and Mendell 
( |2000|) have shown that it does not damp out the r-modes (except for a very small subset of 
the entrainment models employed). But the main point for the present discussion is whether 
or not superfluidity leads to a new set of r-modes, thus making the problem richer than 
the normal fluid one. Unfortunately, this is hard to discern from the discussion of Lindblom 
and Mendell ( |2000| ). While they argue for the existence of some type of superfluid r-modes 
from their numerical results, in the sense that they see the counter-motion of the neutrons 
with respect to the protons, they find that these modes differ little from the ordinary fluid 
r-modes at both the lowest- and second-order in the angular velocity of the background con- 
figuration. Also, they argue that any analogy with the spheroidal type of superfluid modes 
discussed earlier is only superficial, since the fluid motion is not dominated by a deviation 
from chemical equilibrium. We will argue below that a distinct set of superfluid r-modes 
does, in fact, exist. Moreover, we will show that the frequencies of these modes, even at the 
lowest-order in angular velocity, differ from those of the ordinary fluid r-modes, and depend 
in an essential way on the strength of the entrainment effect. These results shed light on 
some of the (so far unexplained) peculiarities on the results regarding the effect of mutual 
friction on the unstable r-modes. 

Carter, Langlois and collaborators ( |Carter 1989| ; |Carter fc Langlois 1995a] ; |Carter & 
Langlois 1995b| ; |Uarter fc Langlois 1998 ; Langlois et al 199"8| ; Comer fc Langlois 1994 ) have 



been developing a system of equations that can describe superfluids in general relativistic 
neutron stars. They are the relativistic analogs of the equations used by Lindblom and 
Mendell and Lee for the Newtonian regime. Until recently, the relativistic equations had 
not been used for modeling various scenarios but this situation is rapidly changing. Comer, 
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On the dynamics of superfluid neutron star cores 7 
Langlois and Lin Q1999Q have used the relativistic equations to model the quasinormal modes 



of a simple neutron star model, and were successful at extracting the general relativistic 
version of the superfluid modes. Andersson and Comer ( |2001| ) have built models of "slowly 
rotating" general relativistic superfluid neutron stars. A unique feature of their formalism 
is that the neutrons are not a priori forced to corotate with the protons. We will use the 
Newtonian limit of these equations for the main analysis presented below. One reason for 
this is to show how to connect the general relativistic formalism with the Newtonian one. 
Another underlying motivation is to show how to incorporate the current microscopic models 
of entrainment in the general relativistic case. 

The layout of the paper is as follows: In Section II, we review the standard reasoning 
that leads to a demonstration of p- and g-modes in cold (single fluid) neutron stars. This will 
provide a convenient background for the analysis of the following three sections (III, IV, and 
V) of the perturbations of superfluid neutron stars. Aside from the main points discussed 
above, a practical outcome of these calculations will be a simple, approximate formula for 
the superfluid mode frequencies that has not been given previously in the literature. Largely 
just to simplify the analysis we use the so-called Cowling approximation, i.e. the variation 
in the gravitational field will be ignored. In Section VI we unambiguously demonstrate the 
existence of two distinct classes of r-modes in a superfluid star. A brief summary is given in 
Section VII. Finally, in the Appendix we give the steps used to derive the Newtonian limit 
of the general relativistic superfluid field equations, and comment on the relation between 
our formulae and those used by Lindblom and Mendell. 

2 STELLAR PULSATION PRIMER 
2.1 The standard picture: p/g-modes 

Before turning to the superfluid case, it is useful to recall the "standard" theory of stellar 
pulsation ( |Cox 1980| ; |Unno et al 1989| ). For an ordinary fluid star, we need to consider the 
linearized version of the perturbed Euler equation 
d 2 £ 5p 1 

H = ^VP--V5P-V5$ (1) 

at 1 p 2 p 

(where £ is the fluid displacement vector, P denotes the pressure, p is the density and <3? is 
the gravitational potential) as well as the (integrated form of) the continuity equation 

5p + V-(pg)=0. (2) 
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8 N. Andersson and G. L. Comer 

We are using 5 to denote Eulerian perturbations, while A will be used (in this section only) 
to indicate Lagrangian variations. 

It is customary to introduce the adiabatic index of the perturbations I\ via 

or, in terms of the Eulerian variations, 



PV 

5P = — Sp + i 

p 



p 



- Vp - VP 



PV 

ee —ISp + PT^-A) (4) 

which defines the Schwarzschild discriminant A. 

For spherical stars we can now readily rewrite the Euler equation as 

Once the equation is written in this form we can deduce that the fluid motion is affected by 
(neglecting 5$) two restoring forces: the pressure variation and the buoyancy associated with 
A. The latter is relevant whenever the star is stratified, either by entropy or compositional 
variations. Since we are considering neutron stars we can to good approximation assume that 
the temperature is zero. This means that we can neglect any internal entropy gradients. Still, 
as was pointed out by Reisenegger and Goldreich ( |1992| ), we cannot assume that A = since 
any variation of the internal composition will lead to an effective buoyancy force acting on a 
fluid element. For neutron stars, the contribution due to the varying proton fraction is likely 
to be the most important effect. 

We want to study @ and try to infer the nature of the various modes of pulsation that 
the star may have. In doing this we assume that the Cowling approximation is used, i.e. we 
neglect the variation 5$ in the gravitational potential. If we further assume that the fluid 
element remains in pressure equilibrium with its surroundings (in such a way that SP = 0) 
we have 

AP « £■ VP ee pi- g , (6) 

where we have used the standard definition of the local gravitational acceleration g. If we 
also use the continuity equation, 

Ap = 5p + i-Vp=-pV-i, (7) 
we find that the radial component of the Euler equation becomes 
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= gC A = -NY (8) 

which defines the so-called Brunt- Vaisala frequency N. In other words, we have oscillatory 
motion whenever N 2 > 0. The resultant modes of pulsation are known as the g- modes, as 
they are essentially governed by gravity and internal stratification in the star. Whenever 
iV 2 < the perturbation is unstable. 

We can study the nature of the star's modes in more detail in the following way (cf. 
Reisenegger & Goldreich ( |1992| )). Since g is purely radial for a spherical star the horizontal 
component of the Euler equation (|5|) leads to 

£l = ~^V ± 6P , (9) 

where we have assumed that the perturbation has time- dependence e lujt . Combining this 
with (0) we find that 

where we have assumed that the angular dependence of 5P can be represented by a single 
spherical harmonic Yi m (9,(p) (which is natural since the pressure variation is a scalar) such 
that 

1(1 + 1) 

v ± • v ± 5P = -^t^sp . 



We also have the radial component of (|5J): 

In these two equations we can replace Ap by noticing that 



Ap = pk {5P ~ P9C) (12) 
(having used the fact that g = —ge r ), and we get 



\U*e ) - 

r A or PL i 

and 



1(1 + 1) p 



uj 2 r 2 PY\ 



5_P_ 

P 



(13) 



Wr 6P + PJ\ 6P = {U}2 + 9A)C ■ (M) 
Let us now introduce new variables 

r 2 C 



(15) 
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10 N. Andersson and G. L. Comer 

SP = <f)5P , 

where 

g/c 2 s dr 



cxp 



and the sound speed is 
AP PT l 



C 2 S = 



Ap P 

With these definitions our two equations can be written 
dt r.o oi r 2 5P 



dr 



OJ 



puj 2 c 2 s (t) 2 ' 



where we have introduced the so-called Lamb frequency 

r2 i(i + iy s 



and 
dSP 



dr 



to 2 + gA 



(16) 



(17) 



(18) 



(19) 



(20) 



(21) 



Given this we can reduce the problem to the following ordinary differential equation for 



d J poo • c 
dr I 



2 ^ 2 



T 2 



-id? 
dr 



[u 2 - N'< 



ta2 



= . 



(22) 



From this we can draw some very important conclusions. We see that the problem reduces 
to the Sturm-Liouville form both for high and low frequencies. For large oj 2 we get 



dr I r 



dr 



+ [oj 2 - N< 



i A* 



ri2 



. 



(23) 



Then standard Sturm-Liouville theory tells us that there will be an infinite set of modes 
which can be labelled by the number of radial nodes (n) of the various eigenfunctions, and 
for which oj n — > oo as n — > oo. In the opposite limit, when oj 2 is small, the problem becomes 



d 


\ W> 2 d?\ 


\N 2 ' 
— - 1 


pi r (j) 2 


dr 


[/(/ + !) dr j 


ur 









(24) 



and we deduce that there will be another set of modes, with eigenfrequencies such that 
u n — > as n — > oo. 

These sets of modes are the p- and g-modes, respectively, and we can estimate their 
frequencies in the following way: Assume that the perturbations have a characteristic wave- 
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length k~ l (such that the various functions are proportional to exp(iAr)). Then we can 
readily deduce the simple dispersion relation 

^ 2 = ^(iV 2 -c 2 )(LjW). (25) 

Here we must have u 2 > for stability, and we see that we will have pulsating modes (k 2 > 0) 
in two different cases. Either u 2 must be smaller than both N 2 and Lf, or it must be larger 
than both these quantities. For I » fcr we estimate the mode-frequencies as uo 2 pa L 2 for 
the p-modes and u 2 ~ N 2 for the g-modes. 

So far, we have discussed familiar textbook results ( Cox 1980 ; [Unno et al 1989|) . We did 



this in order to be able to compare and constrast these results with the corresponding ones 
in the superfluid case. This is the aim of the rest of the paper. 

2.2 The first step towards two fluids 

As a first step towards considering a superfluid star we approach the above pulsation problem 
in a somewhat indirect way. Let us suppose that the star is composed of two distinct species 
of particles, which we will think of as the superfluid neutrons in the core and the "protons" 
as defined in the Introduction. We denote their respective number densities by n n and n p . 
Then standard thermodynamical considerations tell us that the pressure can be determined 
from 

dP = n n dp n + n p dfi p , (26) 

where pi (i = n,p) are the two chemical potentials (which do not include the rest-masses). 
In general, there should also be a term proportional to d\v n — v p \ 2 , where v n and v p are 
the neutron and proton velocities (cf. the Appendix), respectively, but since the background 
velocities are zero (or at least the same for both fluids) in the cases we will consider, such a 
term will not contribute in what follows. By introducing 

tin = , tip = , (27) 

m n m p 

we get the following relation between the Eulerian pressure variation and the variations in 
the two chemical potentials: 

5P = p n 5fl n + p p Sfl p , (28) 
where p n = m n n n and p p = m p n p . 
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12 N. Andersson and G. L. Comer 

Chemical equilibrium for our system corresponds to \x n = \i p . It is natural to introduce a 
variable that describes the deviation from equilibrium introduced by the fluid motion. Thus 
we define 5/3 by 

5 fi p = 5fi n + 5(3 . (29) 

For simplicity, we have assumed that the two particle masses are equal: m n = m p = m. To 
define 8(3 in this particular way (in terms of the tilde variables) may seem peculiar, but it 
is useful since this variable then has exactly the same meaning as 5/3 in the series of papers 
by Lindblom and Mendell (|1994fc |1995|) . Anyway, the above relations enable us to write 

6fi n = — - ^5(3 . (30) 
P P 

The corresponding thermodynamic condition is that 

da n = -dP - ^d(3 . (31) 
P P 

From this we can immediately read off that 
I ( dpm\ P P _ ( dp r , ' 



P \ dP J(3 P K 9 ! 3 ) P 

and using the equality of mixed partial derivatives, we can deduce the useful identity: 



(32) 



Let us now assume that the motion of each species of particles is determined by the 
variations in the chemical and gravitational potentials according to (see the Appendix for a 
justification) 

^ + V(5fi n + 5$) = (34) 
and 

_£ + V(% + = . (35) 

In view of the earlier relations, the second of the two equations can be written 

^ + V(5jl n + 5(3 + 5®) = . (36) 

Given (0) and (|36|) we can make an important observation: The neutrons and protons 
will only move together (in the sense that 5v n = 5v p = dt£) if 6(3 = 0. Intuitively this is 
obvious, since it simply says that the perturbation keeps the fluid in chemical equilibrium 
if the two species move together. Yet, it will prove a useful observation later. We note 
that if the two particle species move together we retain the standard Euler equation @ 
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in the particular case A = 0. This is the first hint of a result that has repercussions on 
the pulsation properties of superfluid stars. As we will now show, the case of barotropic 
perturbations (A = 0) in the standard description, corresponds to 5/5 = in the two-fluid 
picture. 

We elucidate this correspondence in the following way: Introducing a new variable cor- 
responding to the average velocity 

^± = —5v n + —5v p , (37) 
at p p 



(i.e. p dt£+ is the total mass-density current) and combining fl34|) and (|36|) , we readily get 

-|± + V5p n + V5$ + ^V5(3 = (38) 
eft p 

or 

§ +V (^) + V«-^(a)^.o. (30, 
Comparing (|3~9D to the standard Euler equation fl5|), and identifying £ + = £, we see that we 



must have 

^feV^T 1 ^^' (40) 

In other words, the one-fluid Schwarzschild discriminant is intimately linked to the variation 
5(3 in the two-fluid picture, and indeed 5/3 = implies A = 0. This is quite natural since (in 
the case of a chemical composition gradient) the Schwarzschild discriminant describes the 
extent to which a given perturbation drives a fluid element away from chemical equilibrium. 
This observation implies that the pulsation properties of a two-fluid model with 5(3 = 
should be analogous to the case of barotropic perturbations for which there are no non- 
trivial g-modes. In other words, in this case one expects only to find a set of p-modes. 

3 THE EQUATIONS GOVERNING A PERTURBED SUPERFLUID 

We now consider the full two-fluid model for superfluid neutron stars, i.e. consider the super- 
fluid neutrons as being coupled to a conglomerate of charged particles (protons, electrons, 
nuclei etcetera). Adopting the notation of the discussion in the previous section we first 
rewrite the linearized form of the equations ( |A24j )-( |A27| ) from the Appendix in terms of 



a set of variables that are intimately linked to the physics of the system. The first two of 
these are the average velocity <9 4 £ + and the deviation from beta-equilibrium 5(3 that we in- 
troduced in Section 2. Together with these we use the pressure variation 5P and a variable 
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reflecting the relative motion of the neutrons and the protons. In order to retain complete 
correspondence with the equations used by Lindblom and MendellQ we define 

ir-aS (41) 

where we note that the case = corresponds to the two species of particles moving 
together. We have introduced 

det p = PnnPpp - p 2 np , 

Pn Pnn Pnp > 

Pp = Ppp + Pnp ■ (42) 

The mass-density matrix element p np is what allows for the entrainment effect (|Andreev & 
Bashkin 1975|) , and it can be rewritten as 



Pnp = -ep n , (43) 
where e is the so-called entrainment parameter. We also have 

Pnn = Pn(l + e) , Ppp = Pp + £Pn ■ (44) 



In the model used by Lindblom and Mendell ( |2000|) one has 



e=M--l h (45) 
Pn \m* p J 

where m* is the effective proton mass (which enters because the protons form a Fermi liquid 
and it is thus their associated quasiparticle features that are paramount, see Sjoberg (|1976|) ). 

The first continuity equation (i.e. the sum of the linearized versions of (|A24|) and ( A25|) 
of the Appendix) takes the form 

6p + V (46) 

and we also arrive (by adding the linearized versions of ( |A26|) and (|A27 ) of the Appendix) 
at the same equation of motion as in Section II: 

#-(tB(I)/^°- 

* This is obviously not necessary but we want to avoid unnecessary confusion among readers who are well acquainted with the 
relevant literature. 
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We have also assumed (for reasons of clarity) that the Cowling approximation is made, i.e. 
5$ = 0. 

In addition (by subtracting the linearized versions of ( |A26| ) and (|A27 ) of the Appendix) 
we have a second equation of motion relating £_ and 8(5: 



dt 2 



+ V5f3 = 



(48) 



and (after some manipulations^) the final continuity equation can be written 

p ,d( Pj \ SP + ^d (pA 6p+ 



OP 



P 



P 



P 



pL ■ v f & 



V 



p d/3 \p n/ p 
( det p 



P 







(49) 



These four equations are identical to those derived and used by Lindblom and Mendell ( |1994| ; 
19951) . 



In these equations it is, however, difficult to discern the meaning, and thus importance, 
of the various thermodynamic derivatives that appear. Fortunately, the situation can be 
clarified somewhat by introducing the two local sound speeds that our two-fluid system of 
equations admit (just like their mathematical twins, the superfluid equations constructed for 
superfluid helium by Landau ( Putterman 1974]) ). Using the standard technique of expanding 
the fluid variables as plane waves, and considering the background fluid to be homogeneous 
and at rest, we find that these sound speeds are obtained as solutions for u 2 from the 
following quadratic: 







u 



m 



u 



P„2 



+ 



Pp 

Pn 



2c 



Tip 



2 , 2 2 
U +C n C p 



1 + Pp 

Pn 



m r . 



ml 



-'np 



~"np 



(50) 



where 



Pn 



dfin 
9p n 



t The derivation of this equation is somewhat involved. The first step is to "invert" the velocities to find 

8t V PPn 

Sv p = » U+ + — £- 

at \ pp P 



Next one inserts the relationships obtained from (A37) into the linearized form of the proton mass-continuity equation, using 
also the inverted velocities given above. Finally, one uses the linearized form of the neutron mass-continuity equation to obtain 
an equation for 9j&. 
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Pp 



dp P 



Pn 



(51) 



OPp dp n 

When the ratio p v j p n is small we see that the sound speeds are essentially c 2 and (m p /m*)Cp. 
In terms of these definitions we can now write the various thermodynamic derivatives in the 
(still exact) forms 



1 

IT 



X y 



X 



I x p 

X f, 



X 



Cnp 



\ _|_ ZJL I — I 2—^- I " np 



^np 



X 



X ■;■ 



X ., , 



(52) 



' dp^ 



X r 



1_ Z 



Xp 
2 



■up 



^np 



X 



d 



p P 



X ;■ 



dp \p r , 



I X P 



1 + 



X , 



X .,■ 



I — ZJL ) I " np 



^np 
(-a 

2 



(53) 



-1 



X 



+ 2- 



. Xr 



X 7 



u 1ip 



and 



A 



c 2 c 2 

P X n Xp 



X; 



2 "F | ~ n P 1 | ~ ,lj l J 

X ■;■ 



(54) 



(55) 



where x niP = p n ^j p and we have used ( |A32|) - (|A35|) given in the Appendix. (Also, one should 
not confuse the use of the A symbol here with its earlier role as a Lagrangian variation.) 
The utility of these expressions is that they will illuminate the role of the proton fraction in 
determining the order of magnitude contributions from individual terms in the dispersion 
relations that will be written below. 

^)-(|55]) above, and which will be used later, is 

d 



An identity that follows from 

2 



eq 



dp 



Pi. 



Pp 



(56) 



In anticipation of later results, and in order to facilitate comparisons between the ordinary 
fluid and superfluid cases, we define for the superfluid system of equations two frequencies: 



© 0000 RAS, MNRAS 000, 000-000 



On the dynamics of superfluid neutron star cores 17 



1(1 + l)cL 

$ ee r / eq (57) 
and 

The first of these corresponds to the standard Lamb frequency, while the second has (as will 
be discussed later) similar character to the Brunt- Vaisala frequency. 



4 OSCILLATIONS OF NONROTATING SUPERFLUID STARS 

We now want to analyze the superfluid perturbation equations in a way similar to that used 
in Section II for the standard pulsation problem. There are several reasons for doing this. 
The most obvious one being that we want to contrast the pulsation properties of a superfluid 
star with the standard ordinary fluid results in order to see whether future observations may 
be able to distinguish between the two cases. From the general nature of the equations one 
might expect that the character of the various modes of oscillation may be rather different in 
the two cases. After all, in the superfluid case we are allowing the two fluids to move more or 
less independently and so would seem to have brought in additional fluid dynamical degrees 
of freedom. On the other hand, the g-modes in the standard case depend crucially on the 
stable stratification mainly associated with the varying chemical composition ( [Reisenegger 
fe Goldrcich 199^ ). It is thus not at all clear what will happen to these modes if we allow 
the neutrons and the protons to move relative to one another. 

4.1 Two simple limiting cases 

Just as in the ordinary fluid case it is interesting, and potentially instructive, to consider 
what happens if we freeze various degrees of freedom. For example, the equations governing 
£ + and 5P reduce to the standard one-fluid equations (cf. Section II) and it seems reasonable 
to think that these variables should therefore reflect the ordinary fluid properties. At the 
same time, we have seen that 5/3 was to a certain extent accounted for in the standard 
picture via the buoyancy and the Schwarzschild discriminant A. Finally, the variable is 
clearly unique to the two-fluid system and could therefore be expected to bring some new 
features to the pulsation problem. 

Let us first consider the case when the neutrons and protons move in such a way that 
they remain in chemical equilibrium. This would correspond to 8j3 = 0, and we immediately 
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see that this requires = as well. In other words, the neutrons and the protons must 
move together (quite intuitively). From the analysis in Section II (cf. the discussion that 
leads to (0)) we already know that the two equations (J46[) and ([47]) reduce to the standard 
equations for barotropic (nonstratified) stars (A = 0). Hence, the only non-trivial modes we 
expect to find are the p-modes. In this case one can easily show that they will have frequency 

ujI » A 2 (59) 

noting, however, that c 2 is being replaced by cL, i.e. the sound speed defined in (|52|) . This 
is, of course, exactly the ordinary fluid result of Section II in the limit A — ► 0. 

As a side remark it is interesting to note that one would not normally expect the above 
assumptions to lead to the two remaining equations, flUf) and (f49l) , also being satisfied. 
However, in the present case (|48"D is trivial and, for particular models, one can also satisfy 
((49p. This requires 



(60) 



A(V| =0 . 

In other words, we need to have p p oc p or p p = 0. The latter does, of course, correspond to 
a star in which there is only one species of particles so the equations of motion must reduce 
to the standard ones for a single fluid. 

In Section II we saw that the nature of the g-modes could be deduced by requiring 
that the motion be such that the fluid remains in pressure equilibrium, 5P = 0. Let us 
consider the equations from the previous section in this case. Conveniently "forgetting" (for 
the moment) the conservation of mass equation ( |46"D we need to consider 
d5(3 



uj 2 p 



dr r 2 det p 
and 







(61) 



dZ_ 

dr 



'detp 2 
c 4 A 

r 2 c 2 A 
— ^— 5f3 
P 





—Af 


')] 




_ 9 





X 



(62) 

which follows from (|47|)-(f49D, once we separate the radial and horizontal components of the 
equations (as in Section II) and introduce the new variable (£H is the radial component of 

f-) 

Z = ^l r 2 C , (63) 



as well as use the various definitions and relations given at the end of Section III. We now 
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combine the two equations to find 

1 2 N 



d 
dr 



P 



r 2 c 2 A 

eq 



1 + 



—M 

9 



'detp 2 



A" 2 ] ' , 



dr 



uj 2 p 



, ] ,- Z = (64) 
r z det p 

and it is easy to see that, just as in Section II, we have a Sturm-Liouville problem for large 
uj 2 . Hence, we expect there to exist a set of modes for which u 2 — > oo for the large overtones 
(with index n — > oo). 

To elucidate the nature of these modes, let us derive the appropriate dispersion relation. 
Assuming that the perturbation variables behave as exp(ikr) we readily arrive at 



k 



c 2 A 
det p 







to 


[('♦ 








_ 9 





detp C 2 + Af 2 



c 4 A 

v ^eq 



and in the limit when I » kr we find 



uj 2 ~ 

s 











—M 






_ 9 





detp C 2 + Af 2S 



c 4 A 



(65) 



(66) 



Assuming that the sound speeds are well-behaved if x p is small, then we can infer from (|52|) 
that Cg q , and hence Cf, are roughly independent of the proton fraction to leading order. We 
also infer from ( ^2[ ) and ( |55| ) that the combination detp/A is independent of the proton 
fraction to leading order. Finally, we see from (|55| ) and ([5^) that 

i (dp\ 2 

which implies that 



(67) 



N" 2 



X 



p 



Thus, to linear order in the proton fraction we see 

1 



10 ~ 



$ 1 



—M 

9 



det p 
c 4 A 



+ Af 2 . 



To the lowest-order in the proton fraction, the dominant contribution to uj s is 

A 2 



2 _ .2 det P _ + 1) 2 



C 4 A 



r z 



(68) 



(69) 



(70) 



Hence we deduce that these "superfluid" pulsation modes are essentially governed by the 
proton sound speed c 2 , and also that the entrainment effect plays an important role in 
determining their frequencies. We also notice that contributions due to M are negligible at 
this order. This approximate result provides a natural explanation for the result of Comer 



© 0000 RAS, MNRAS 000, 000-000 



20 N. Andersson and G. L. Comer 

et al ( |1999| ) that the ordinary and superfluid frequencies are virtually the same when the 
adiabatic index of the neutrons (treated as a relativistic polytrope) is set equal to the index 
for the protons (also treated as a relativistic polytrope). In that case, the sound speeds of 
the two fluids are virtually the same. Furthermore, in the models considered in Comer et al 
there is no entrainment since m* is set equal to m p . 



4.2 The fully coupled case 

Having established the existence of two distinct classes of pulsation modes in a superfluid 
star we now return to the full problem, and consider the modes of pulsation for the system 
of equations ©-(H). 

Proceeding in the now familiar way, we separate the radial and horizontal parts of (f47j) 
and fl48|), and then combine the results with the two conservation laws ([46]) and (f4~9"D. This 
way we arrive at a set of four equations 



r 2 dr ( Cg q ^ + ~'~ 
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IT 



1(1 + 1) 

uj 2 r 2 



5P 

— + 

p 



i (dp; 
p 



5(3 = 



+ ^5P - P ur + + 9 
dr c 2 eq 



5(3 = 



r 2 det p d 
p dr 



5(3 = 



1 dZ 

r 2 dr 



p\ d (p p \ 1(1 + 1) det p 

p 8(3 \p n ) 



ijj 2 r 2 



P 



5(3+ 



5P 

gt 

p 



o 



(71) 



dp) 

where is the radial component of £ + and the variable Z was defined in Section IVA. 

As in Section II the first two equations simplify considerably if we introduce the inte- 
grating factor as defined by ( |T7| ) (replacing c s with c eq , of course), and then work with the 
variables 5P = <p5P and Y = r 2 £^/0. If we also introduce the characteristic wavelength of 
the pulsation through exp(ifcr), our first two equations take the form 



ikY + 
and 



1(1 + 1) 



-eq 



5P 



r 
P 



'dp\ 5p 



(72) 
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(73) 



These two relations can easily be combined to provide a relation between 5P and S /3: 
l{l + l) + k 2 r 2 2 



5P 



[ikg + u 2 ] x 
2 .( d P 



-eqi 



(74) 



After some similar manipulations, one can show that the two remaining perturbation equa- 
tions lead to 



[lj 2 - ikg] 



' dp" 



5P 



1(1 + 1) + k 2 r 2 



det p— 



2 2 9 (Pp" 



(75) 



These two equations can be combined to provide the dispersion relation for waves propagat- 
ing in a superfluid star 



[- 4 + *Vk 2 q (f^ 



Kl + 1) + k 2 r 2 2 



" C eq - W 



det p - w p n 



<9/3 \Pn / 



+ 



' dp s 



(76) 



In principle, this equation contains information equivalent to the dispersion relation of Lind- 
blom and Mendell (PM|) (cf. their equation (82)). 

We can do one final rewriting of ( f76|) as a quadratic in tu 2 , which we could then solve. 
Using the various definitions and relations listed at the end of Section III, the fully coupled 
dispersion relation takes the final form 



= uo A - 



k 2 r 2 



1 + 
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1(1 + 1)_ 

k 2 r 2 \ 
1(1 + 1)) 

k 2 r 2 



1 + 



det p 
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A 2 + 



—A 

K 9 , 



Af 2 uo 2 + 



det p 



c 4 A 



C 2 +N 2 



(77) 



It is, however, not particularly instructive to write down the formal solution to this equation 
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It is much better to first simplify it somewhat. This provides a better insight into the relevant 
physics of the solution. 



4.3 Results for a small proton fraction 

Our aim now is to use the results of the previous section to infer how the presence of a 
superfluid affects the pulsation modes of a nonrotating neutron star core. In principle, we 
have drawn the main conclusions already. There will exist two more or less distinct classes of 
modes. One of these comprise the standard p-modes, while the other has unique properties 
(although see comments below) due to the presence of the superfluid. Both these families 
of modes are such that the eigenfrequencies -> oo as the mode- number n — > oo. In 
other words, the two sets of modes will be interlaced in the pulsation spectrum of the star. 
These are, however, only qualitative results and it would clearly be useful to make more 
quantitative estimates. 

The key point to simplifying the fully coupled dispersion relation is that M 2 ~ x p which 
means that the frequencies can easily be determined to linear order in the proton fraction. 
Doing this (and considering also the limit where I » kr) yields the two solutions 



$ 1 + 



1 - 



det p 
c 4 A 

cq . 



-1 



c, 
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detp 2 
^A ' 
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1 det p 


-1 
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c 4 A 




_ 9 





+ Af 2 



(78) 



(79) 



To lowest-order in the proton fraction, we recover the ordinary, i.e. u 2 , and superfluid, i.e. 
uj 2 , solutions respectively. 

Let us now address one of the main questions that motivated the present investigation: 
What happens to the g-modes when the neutron star core becomes superfluid? We recall 
the discussion by Reisenegger and Goldreich ( |1992j ) that showed that a non-superfluid neu- 
tron star will have a distinct set of g-modes whose nature is determined by the chemical 
composition gradient associated with the varying proton fraction. From Section II we know 
that these modes will have eigenfrequencies such that u 2 — > as n — > oo. We have already 
shown that there will not be a family of modes with this character in the superfluid case. 
Thus, even though we clearly have a varying proton fraction, a superfluid neutron star core 
will have no pulsating g-modes, in accordance with the numerical results of Lee ( |1995| ). This 
is part of the answer to the question posed above, and in order to provide the rest of the 
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answer, we will need to consider the space of time-independent solutions to the linearized 
equations, i.e. the zero-frequency subspace. But before doing that, we want to emphasize 
some consequences of having a varying proton fraction. 

The effect of a composition gradient will be analyzed as follows: In the ordinary fluid case 
one can estimate the Schwarzschild discriminant due to a variation in the proton fraction as 
(see Unno et al (|1989|) and the discussion in Lee (|1995|) that surrounds Lee's equation (47)) 

A =-\7T-\ er-Vp p (80) 



P \dp P , 



p 



(assuming the star has zero temperature). This then leads to g-mode frequencies 
u 2 nN 2 = -gA=- 9 -(^-) e r -V Pp , (81) 



or 



i2 (dp_\ (dpA (g_ sl 



^-mii^Mtj • <82) 

Now consider the dominant contribution to the superfluid mode frequency (|79|) when x p and 
c 2 are both small; it is simply 

■ 2 -M 2 ^xJ^\ . (83) 



s 



\ P 



From this we see that uj 2 ~ N 2 « M 2 ~ uj 2 s when the proton fraction and c 2 are small. Thus 
we conclude that, even though a superfluid neutron star core does not support a distinct 
class of pulsating g-modes, the buoyancy associated with internal composition gradients can 
be relevant. Of course, the effect depends crucially on the proton fraction and cannot be 
distinguished unless c 2 is small. 



5 THE ZERO-FREQUENCY SUBSPACE 

We will now focus our attention on the perturbations of the superfluid equations that belong 
to the zero-frequency subspace. That is, after perturbing (and linearizing) the equations, we 
assume that the velocity perturbations, as well as Sn n and Sn p , are time- independent (cf. 
( Lockitch fc Friedman 1999| ; |Lockitch et al 2001|) ) . We then get from (|A24l) and ( |A25|) : 

= d l (pjv l n ) , = d i ( Pp 6v p ) . (84) 

Meanwhile, the two linearized Euler equations and the equation for the linearized Newtonian 
potential lead to 

= d i (5$ + 6ji n ) , = d i (6$ + 6ji p ) , (85) 
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and 

V 2 5$ = ArtG (5p n + 5p p ) = AnGdp . (86) 

From these equations we immediately deduce that we can split the perturbations into two 
independent sets. The first set has the perturbed velocities vanishing but the other per- 
turbations non-zero. The second set has 8v l n and 5v l p nonvanishing while the local matter 
perturbations 5p n = 5p p = Sjl n = 5jl p = 5$ = 0. The general solution will be a superposition 
of solutions from each set. 

5.1 Neighbouring spherical equilibria 

We begin by interpreting the meaning of the first set of perturbations distinguished above. 
By taking the difference between the two equations in (|85|) it is easy to show that we must 
have 

diS/3 = . (87) 

Given this, we can also show that the sum of the two equations leads to 

diSP = -5pdi<$> - pdiS® . (88) 

These two equations follow simply from the perturbed versions of the equations that deter- 
mine the unperturbed configuration, i.e. 

d l ($ + p n ) = , d l (®+p p )=0, (89) 

using also the earlier equation that relates pressure variations to chemical potential varia- 
tions. And just as for the equilibrium equations we can show that the solutions must have 
spherical symmetry. Thus, the interpretation of this set of perturbations becomes clear: They 
simply correspond to a neighbouring spherical equilibrium model. Provided that 5(3 = the 
perturbed model is in chemical equilibrium. 

5.2 Non-trivial velocity perturbations 

We now turn to the other set of perturbations, for which the three-velocities are non-zero 
and represent convective currents. We first have to account for the fact that perturbations 
of a spherical star can be decoupled into two different classes: spheroidal and toroidal. 
Spheroidal velocity perturbations take the form 

5v 3 = - (Wj(r), Vj(r)d e , Y lm , j=n,p (90) 



r V sin 6* 



© 0000 RAS, MNRAS 000, 000-000 



On the dynamics of superfluid neutron star cores 25 

where we have used a standard decomposition in terms of spherical harmonics. With these 
definitions we readily get from (|84j) 



= 1(1 + l)n n V n - — (rn n W n ) , (91) 

= 1(1 + l)n p V p - ^ (rn p W p ) . (92) 

Foregoing here a detailed discussion of boundary conditions, we can establish that non-trivial 
solutions exist. In particular, we see that there are two functions (W n and W p , say) that can 
be freely specified. Because they are spheroidal, and have vanishing pressure and chemical 
potential perturbations, we identify these solutions as the g-modes "missing" from the time- 
dependent perturbations. In this respect, we see that superfluid g-modes behave like those 
in the barotropic ordinary fluid case, i.e. when the equation of state of the perturbations 
is the same as that of the background ( |Lockitch fc Friedman 1999| ; [Lockitch et al 200 1|) . 



However, the doubling of the modes (because there are two free functions) is unique to the 
superfluid. 

Toroidal three-velocity perturbations can be written as 

80, = - ( 0, -Uj^de) Y lm , j —n,p (93) 

r V sin 6* / 



For velocity perturbations of this form Eq. (p4|) is automatically satisfied for arbitrary U n 
and U p . Since there are two functions that can be freely specified we deduce that there will 
be two sets of toroidal modes in the zero-frequency subspace. When the star is rotating this 
should lead to the presence of two classes of r-modes. We will verify this in the next section. 

6 ROTATING SUPERFLUID STARS: QUALITATIVE INSIGHTS INTO 
THE R-MODES 

In the last two years various aspects of rotating neutron stars have been under intense 
scrutiny following the discovery ( [Andersson 1998| ; |Friedman fc Morsink 1998| ) that the emis- 



sion of gravitational waves drives the so-called r-modes unstable. Initial studies of the prob- 
lem indicated that the r-mode instability might cause a newly born neutron star to spin 
down to a rotation rate comparable to that inferred from observation for the Crab pulsar 
( [Lindblom et al 1998 ; Andersson et al 1999 ), in the process radiating gravitational waves 



that may well prove detectable with the generation of large-scale interferometers due to 
come online in the next few years ( Pwen et al 1998 ). One issue of utmost importance for 
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the r-mode problem concerns the role of superfluidity ( [Andersson fc Kokkotas 200 1| ). So 



for example was it originally thought that dissipation due to the superffuid mutual friction 
would counteract the instability in a significant way. The only available calculation of this 
effect, due to Lindblom and Mendell ( |2000| ), suggests a rather different picture, however. 
Their results indicate that the mutual friction will typically not be able to suppress the 
r-mode instability. But the results also suggest that the effect becomes dominant for certain 
values of the entrainment parameter e. Clearly, this problem is far from well understood at 
the present time. 

Our aim in this section is to study the r-modes at a qualitative level (comparable to our 
study of spheroidal p- and g-modes in the previous sections). This would seem a natural 
starting point for a discussion of the r-modes of superffuid stars, and as we will see it provides 
insights that may well explain various features observed in the numerical work by Lindblom 
and Mendell (pOOO)) . 

To analyze the r-mode problem we will extend Saio's ( |1982| ) vorticity argument to the 



two-fluid problem. Thus we focus on the two Euler equations. When expressed in terms of 
the variables £ + and these equations can be written 

8$£+ + (u ■ V)d t £ + + d t £ + ■ W + V5<1>+ 

V5/i n + x p V5(3 = (94) 

and 

dU- + (u ■ V)dtg- - d t C • W + V6p+ 

2^(<%£_-V)« = 0. (95) 

PnPp 

Here we have assumed that the protons and neutrons corotate in the unperturbed case. The 
relevant rotation velocity is denoted by u in the above equations. We note that it may be 
desirable to relax the assumption of corotation of the two background fluids in order to 
model a real neutron star (cf. the discussion in Andersson and Comer ( |2001| )). The resultant 
problem is, however, much more complicated than the present one and we will return to it 
in future investigations. 

We first translate these equations into the rotating frame. Then we get 

d 2 t i + + 2VI x dt(+ + V<5$ + V5j2 n + x p V5(3 = (96) 
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det p 

PnPp_ 



(d£- ■ V)« = . (97) 

The next step involves assuming that the mode is horizontal to leading order, taking 
the curl of the above two equations and focusing on the radial component of the resultant 



relations. Consider first (|96|) , which readily yields 



V x dti + + 2VxOx^ 







(98) 



since 



V x (x p VSP) = Vx p x VSp (99) 
has a vanishing radial component for a slowly rotating (still spherical) star. We now use 
VxOxC+] r = {fi(V-e+)-e+(V-fi) + 



(e+-v)n-(n-v)e+} 



« (£+ • V)Q . 

L J r 

After these manipulations we have arrived at 



a i (Vx( + ) + 2(( + .v)d 







(100) 



101) 



Finally, we insert in this equation the assumption that the mode we are interested in is 

purely toroidal, i.e. is of the form . 
T 



r 



. -d^-Td^e^Yi 
smv 



l in 



Given this assumption we immediately find that these modes must have frequency 
2mVL 



to, 



(102) 



(103) 



1(1 + 1) ' 

in the rotating frame. This is, of course, the standard r-mode result (|Andersson fc Kokkotas 
2001|) . 

Let us now go through the same exercise for the second Euler equation, (|97|), that appears 
in the superfluid case. The radial component of the curl of this equation can be written (in 
the rotating frame) 

d t {(V x dtC) + 2V x (Q x £_)- 
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2V x 



det p s 

PnP P/ 



. 



(104) 



Here the first two terms are identical to those we encountered in the analysis of fl96|), but 
the last term is new. This term can be written 



2V x ( 7 £_ -U) =2 (V7) x (fi x £_) + 2 7 V x (O x 



where we have introduced 
det p 



7 = 1- 



PnPp 



(105) 



(106) 



It is clear that the first term of the right-hand side of ( |105| ) does not contribute to the radial 
component. The second term, however, does have a non- vanishing radial component. We 
readily find that 

27ft 



[2Vx (jg--H)] r 



r sin 9 



(107) 



d e (sin 9 cos 9£ e _ ) + d v (cos 6C ) 
and if we assume that the vector £_ is toroidal, i.e. takes the form ( |102j ), we have 

[2VX hi--u)] r = ^Jx 
v ' S1I1& 



[d e (cos Od^Yim) - d v (cos 0d g Yi m )] 

= -2imVt-iTY lm . (108) 

We can now combine this result with the final result obtained from ( |96"D to deduce that we 
will have modes with frequency 

^ = Jrn^detp (iog) 



/(/ + 1) p n p p 

From the definition of det p it follows that this is identical to 



2mQ 
1(1 + 1) 



V Pv)\ 



(110) 



where e is the entrainment coefficient. This is a very interesting result since it demonstrates 
the existence of a distinct class of superfluid r-modes, the properties of which are to a large 
extent determined by the entrainment parameter. 

Based on the above analysis we can now discuss the general nature of the toroidal modes 
of rotating superfluid neutron star cores (even though we should advice some caution since 
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( P-10| ) is r-dependent and therefore strictly speaking only describes an r-mode in a thin shell). 
We have seen that, just like in the non-rotating case, there will be two distinct classes of 
modes. The first corresponds to fluid motion such that the neutrons and the protons flow 
together, and the resultant modes are analogous to the standard perfect fluid r-modes. These 
are the modes that the calculations of Lindblom and Mendell ( |200(J| ) were aimed at studying. 
For the second class, the neutrons and the protons are counter-moving. This class of r-modes 
has not previously been discussed in the literature (although see the discussion below). It 
is interesting to note that the two classes of modes are (to leading order) degenerate in 
the absence of entrainment. This is, of course, rather intuitive since the two fluids are then 
effectively uncoupled and the Coriolis force affects each fluid separately. 

We note that the existence of two distinct classes of r-modes in superfluid stars provides 
a likely explanation for the somewhat peculiar results obtained by Lindblom and Mendell 
(|2000|) regarding the effect of mutual friction on the r-mode instability. In their study Lind- 
blom and Mendell found that the mutual friction dissipation timescale was sufficiently long 
that this effect would not suppress the unstable r-modes for most values of e. But they also 
discovered that mutual friction would become dominant for certain values of e, cf. their 
figure 6. Furthermore, they noted that while the variable 5(3 was small in most cases, it 
became large for the particular values of e at which mutual friction was found to be im- 
portant. Given our present conclusions regarding the r-modes in superfluid stars, we can 
interpret the Lindblom-Mendell result in the following way: There are two distinct classes 
of r-modes in a superfluid star. One of these classes, the one for which the neutrons and 
protons flow together, is typically not rapidly damped by mutual friction. The other class of 
modes, however, is such that the two fluids are counter-moving. Since mutual friction tends 
to damp relative motion between the two fluids these modes will be strongly affected by 
mutual friction. This then makes the Lindblom-Mendell results quite natural. For certain 
values of the entrainment parameter they have simply found the second class of modes rather 
than the first. 

The situation may, however, be yet more complicated. It is likely relevant to make the 
following observation: In many situations where the mode-spectrum is studied for a family 
of stellar models for which a single parameter is varied (such as the crustal shear modulus) 
one finds that the modes undergo what are known as "avoided crossings." So for example 
have recent studies of neutron stars with an elastic crust shown that the r-modes associated 
with the core fluid show avoided crossings with the toroidal shear modes in the crust for 
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certain values of the rotation rate. Avoided crossings between pulsation modes are, in fact, a 
regular occurence in studies of more complicated stellar models (|Cox 1980| ; |Unno et al 1989| ) . 
It is interesting to note that, well away from the crossing point, the different modes have a 
distinct nature, but as one approaches the avoided crossings the two modes become similar 
in nature and when one reaches "the other side" they have completely exchanged properties 
(the core r-mode has become a crust shear mode and vice versa). In view of our qualitative 
analysis and the numerical results of Lindblom and Mendell, it would seem natural to predict 
that a more detailed investigation of the superffuid r-mode problem will unveil analogous 
avoided crossings between the co- and counter-moving r-modes as the entrainment coefficient 
is varied. 



7 CONCLUSIONS 

In this paper we have compared and contrasted the pulsation properties of a superffuid 
neutron star core (represented by a relatively simple two-fluid model) with the familiar 
(textbook) results for a normal fluid star. This study provides a qualitative understanding 
of the nature of the oscillation modes of a superfluid neutron star, and thus provides a 
theoretical fundament for various numerical results in this area ( [Lindblom fc Mendell 199~T| ; 



Lee 1995) ). We have shown that a non-rotating superfluid neutron star core exhibits two 
distinct families of pulsation modes. Our approximate results provide a deep insight into 
the physics of these modes, with one set of modes corresponding to the two fluids moving in 
"lock-step" while the second family (the modes that are unique to a superfluid) correspond 
to the two fluids being counter-moving. We have shown that both sets of modes are crucially 
governed by the acoustic properties of the two fluids. This is in clear contrast to the normal 
fluid case, where the two families of modes, the p- and g-modes, are acoustic and governed by 
buoyancy (mainly due to chemical composition gradients), respectively. We have investigated 
what happens to the g-modes as the star becomes superfluid and have confirmed Lee's 
numerical results that there are no propagating g-modes in a superfluid core. In addition, 
we have shown that the "missing" g-modes can be found in the zero-frequency subspace. 
Finally, we have shown how the various modes are affected by the parameters of entrainment. 

Basically, our results are important for two reasons. Firstly, they provide a clear and 
concise description of the problem and the general nature of pulsating superfluid stars. Thus 
they fill what we perceive as a gap in the existing literature. In particular, we believe we have 
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resolved some open questions regarding g-modes in the superfluid case. Secondly, the results 
indicate that future observations of neutron star mode oscillations, eg. by a highly sensitive 
generation of gravitational wave detectors, could potentially help constrain the parameters 
of the large scale superfluidity that is believed to exist in the core of mature neutron stars. 
We will discuss this exciting possibility elsewhere. 

In addition to studying non-rotating stellar cores, we have presented a qualitative analysis 
of the r-modes in a slowly rotating superfluid core. We have demonstrated that there will 
be two distinct families of such modes, and argued that the interplay between these modes 
(eg. via so-called avoided crossings as the entrainment parameter is varied) may provide an 
understanding of puzzling numerical results obtained by Lindblom and Mendell (|2000|) . Our 
current analysis clearly shows that the r-mode problem for superfluid stars is likely to be 
much richer than has previously been appreciated. In particular since, in addition to the 
r-modes, we anticipate that there will be a large class of "hybrid" inertial modes ( [Lockitch] 
fe Friedman 1999] ; |Lockitch et al 2001]) . No attempts to study such modes in a superfluid 
star have yet been made. A better understanding of these various issues is clearly needed 
in view of the fact that the r-modes have been shown to be unstable due to the emission of 
gravitational waves, and that this instability may govern the spin-evolution of neutron stars 
during various phases of their lives. We therefore plan to investigate this problem in greater 
detail in the near future. 
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APPENDIX A: NEWTONIAN LIMIT OF THE GENERAL RELATIVISTIC 
SUPERFLUID EQUATIONS 

The purpose of this Appendix is to derive the Newtonian version of the general relativis- 
tic superfluid field equations. There are several motivations for doing this. For conceptual 
reasons we want to be able to compare our relativistic calculations ( |Comer et al 1999| ; [An-| 
dersson et al 2001| ) to previous results in this research area (in particular those of Lindblom 



and Mendell (|1994| ; |1995| ; |2000|) ) most of which were obtained in the Newtonian regime. We 



also feel that it is useful to address various problems on a qualitative level before working 
out the detailed answers. Given that Newtonian calculations are often simpler, it is natural 
to address the Newtonian problem first. In addition to this, we hope that a study of the 
Newtonian limit will provide insight into ways of extending various Newtonian entrainment 
models to the relativistic regime. 

Al The general relativistic formalism and Newtonian limit 

We begin by recounting the formalism that has been used to model general relativistic 
superfluid neutron stars ( |Carter fc Langlois 1995a| ; |Carter fc Langlois 1995b| ; [Langlois et al 



1998 ; Comer et al 1999] ; [Andcrsson fc Comer 200 1| ). The fluid dynamical degrees of freedom 



are described by n%, the conserved neutron number density current, and n£, the conserved 
proton number density current. The fundamental scalar in the formalism is the so-called 
"master" function A, which is a function of the three scalars r? n = —n n p n p w n 2 = —n v p n p vl and 
x 2 = —n^n?, where n™ = g^n^ and n p p = g^n^. The quantity —A corresponds to the total 
thermodynamic energy density of the entire fluid. 

There are two momenta /i M and which are dynamically, and thermodynamically, 
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conjugate to and n£. They are defined by a general variation (that keeps the spacetime 
metric fixed) of A(n 2 , n 2 , x 2 ), i.e. 

5A = fi p 5n p n + Xp Sn p p , (Af ) 

where 

zv^x + ^X » xm = ^X + /c "X' ( A2 ) 

and 

dA dA dA 

K np = -jr- 2 , K, nn = -2?— , K*> = -2% li . (A3) 

ox 1 an l n on l v 

In this general relativistic context the entrainment effect is seen in that the momentum /i^, 
say, is a linear combination of n p and n£. The generalized pressure \I/ is given by 

* = A - n> p - nP Xp . (A4) 

The equations of motion consist of two conservation equations, 

VX = , V„n£ = 0, (A5) 

and two Euler type equations, which can be conveniently written in the compact form 

<V[^„] = , n^V^] = 0, (A6) 

where the square braces '[ ]' indicate antisymmetrization on the indices. 

In order to derive the Newtonian limit, the general relativistic field equations will be 
written to order c° where c is the speed of light. Formally, the Newtonian equations will 
then be obtained in the limit that the speed of light c becomes infinite. The gravitational 
potential, denoted $, is assumed to be small in the sense that 

- 1 « % < . (A7) 
To order c° the metric can be written as 

ds 2 = -c 2 (l + ^p) dt 2 + 5 ij dx i dx j , (A8) 

where the x % {% = 1,2,3) are Cartesian-like coordinates. 

In the superfluid field equations there are two four-velocities that must be considered. 
This means that we must take into account the fact that the fluids define two different 
proper times: one for the neutrons, to be denoted r n , and the second for the protons, to be 
denoted r p . The two different fluid trajectories are then obtained from the functions 

<(r n ) = (t(r„),<(r n )) , 
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= (t(r p ),x;(r p )) . (A9) 
The respective four- velocities of the neutron and proton fluids are thus given by 

< = f • -<-g- < A1 °> 

Because of the choice of coordinates, the four- velocities satisfy u'^u^ = — c 2 and u^u^ = — c 2 , 
where w™ = g^u v n and = g^u p . 

We will turn this around and use a global time t(r n ) = t(r p ) = t as the parameter for 
both curves. In this case, the two proper times will be given by 



dsl = c 2 ^ = c 2 (l + ^-^^]A 2 , 



-dsl = c 2 rfr 2 = e 2 (l + 2* - & , (All) 

where v l n = dx l n /dt and t>* = dx p /dt are the Newtonian three- velocities of the neutron and 
proton fluids, respectively. Each three-velocity is considered to be small in the sense that 



lOnl « i 



« 1 . (A12) 



c c 

Hence, to the correct order the four-velocity components are given by 
and 



n= 1 -^ + 2^»< = «n, (A13) 



$ V 2 



«t =i -^+2^' «i=«i' ( Ai4 ) 



where f 2 = 5ijV l n v j n and w 2 = 5ijV p v J p . 



Note that the two particle number currents are now written as 
< = n n «/c) , n£ = np(u£/c) . (A15) 
To the correct order, one finds for x 2 = —n^n^ that 



w 2 



x 2 = n n n p \1 + — j , (A16) 
where 

w 2 = 5ij (< - v p ) (vi - vi) . (A17) 

In order to write the Euler equations with the terms to the required order, it is necessary 
to explicitly break up the "master" function into its mass part and internal energy part E, 
i.e., to write A as 
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A = - (m n n n + m p n p ) c 2 - E(n 2 w n 2 p , x 2 ) . 

In this context, E is small in the sense that, for instance, 
E 



< 



« 1 . 



(A18) 



(A19) 



It is also convenient to use a different choice for the independent variables that more 
closely agrees with what is used for Newtonian superfluids, which is the triplet of variables 
(n 2 n) n 2 p) w 2 ). Thus, from now on we assume that E = E(n 2 n) n 2 pl w 2 ). 

With this choice, a variation of A that leaves the metric fixed yields the following: 

dA = —{m n c 2 + n n )dn n — (m p c 2 + ji p )dn p — adw 2 , (A20) 



where 



Hn = 



dE 



dE 



dE 



a = 



(A21) 



dn n dn p dw 2 

Now the JC np , JC nn , and K, pp coefficients of the general relativistic formalism are related to 
the coefficients defined above via 
2c 2 



IC np 



-a , 



Tl n Tl p 



m n c 2 + pL n 2c 2 



n„ 



n 



m p c + n P 2c I if 



n r , 



n 



1 + 2?l°- 



u p ,„ p 

In terms of these variables, the pressure P is seen to be 
P = -E + fi n n n + fipUp . 



(A22) 



(A23) 



The Newtonian limit of the general relativistic superfluid field equations reduce to the 
following set of 8 equations: 



= ^ + d t (npV p ) , 



dt 

= 



v J n dj [ vl 



(A24) 
(A25) 



2a 
2a 



m„n. 



n ll, n 



+ 



m r . 



2a 



s ij tki (v l p - v l n ) a 



-J V n , 



(A26) 
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and 
= 



d ( i 2a 



dt 



2a 



+ 



+s>>9 > + 



2a 



5«5 kl (y l n - v l p ) d jV k p . 



(A27) 



Finally, the equation for the gravitational potential $ is 

<9i<9*$ = AttG (m n n n + m p n p ) . (A28) 
These equations have been derived independently by Prix ([Prix 2001 ), using a Newtonian 



Lagrangian-based variational principle. 



A2 Comparison with the Formalism of Lindblom and Mendell 

In their studies of oscillating superfluid stars Lindblom and Mendell QMendell fc Lindblom 



199TI IMendell 19911; ILindblom & Mendell 19941; ILindblom & Mendell 19951; ILindblom & 



Mendell 2000| ) use the mass densities, p n and p p , and the fluid three- "velocities" V n and V p 
as their main variables. The latter are not the same as the three-velocities used in the main 
text of this paper, rather these "velocities" are the "macroscopically averaged" gradients 
of the phases of the mesoscopic wave functions that describe the superfluid neutrons and 
the superconducting protons and they are actually proportional to the spatial components 
of our neutron and proton momenta. By comparing our mass- conservation equations with 
those of Lindblom and Mendell, which are given by 

dp n 











dt 



+ V- {pnvVn + PnpVp 



' P + V- [p vn Vn + PppVp 



Q t • • \ryn ■ « 1 rvv ■ P j 

where p nn , p pp and p np = p pn were defined earlier in the main text, we find that 



(A29) 



PnV r , 



PnnVn Pnp^p 



PpVp 



P pn^n PppVp 



A final comparison to our fluid momenta yields 



a. 



PnPp 



'Pnp 



2 Pnp PnnPpp 
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This is a very useful result since it enables us to make contact with the model that Lindblom 
and Mendell use to describe the entrainment effect. 

In Section III of the main text, we made reference to equations (69)-(72) of Lindblom 
and Mendell ( |1994| ). Neglecting terms of order (m e /m n ) 2 , where m e is the electron mass, 
then these equations are (recalling that jl n = p n /m n etcetera) 



A 

P 



'dp* 



Pp 



9 Pp. 



+ 



, d P P , 



pn 



(A32) 



'dp) 



A 

P 



' 3p r 



dp n 



Pp 



(Pp ~ Pn) 



Pp 
dpn) 
PPp) 



'dpp 
dpp 



+ 



Pn 



Pn 



(A33) 



d(3 \Pn)p 



A 

P 



2 Ppf9pn) 



d PnJ Pp Pn \ d PpJ p 



f% ( dp P \ 
Pi \ d Pp) Pn 



(A34) 



and 

1_ _ / dp n 
A \dp n 



dpp, 



'dp n 



(A35) 

p p \~ ry / p n \ ~ i- 1> / p n 

Note that these follow by considering p n and p p to be functions of p n and p p , and thus 
re-writing (^) and ( ^9|) in the forms 

5P=(...)5 Pn + (...)5p p , 5(3=(...)5p n + (...)5 Pp , (A36) 
and then inverting to find 

5p n = (...)5P + (...)5(3 , 5 Pp =(...)5P+(...)5(3 . (A37) 
Finally, by introducing the various sounds speeds we arrive at Eqns fl52"D - (|5 
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